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The lack of suitable quantum emitters in silicon and silicon-based materials has prevented the 
realization of room temperature, compact, stable, and integrated sources of single photons in a 
scalable on-chip architecture, so far. Current approaches rely on exploiting the enhanced optical 
nonlinearity of silicon through light confinement or slow-light propagation, and are based on para¬ 
metric processes that typically require substantial input energy and spatial footprint to reach a 
reasonable output yield. Here we propose an alternative all-silicon device that employs a different 
paradigm, namely the interplay between quantum interference and the third-order intrinsic nonlin¬ 
earity in a system of two coupled optical cavities. This unconventional photon blockade allows to 
produce antibunched radiation at extremely low input powers. We demonstrate a reliable proto¬ 
col to operate this mechanism under pulsed optical excitation, as required for device applications, 
thus implementing a true single-photon source. We finally propose a state-of-art implementation 
in a standard silicon-based photonic crystal integrated circuit that outperforms existing parametric 
devices either in input power or footprint area. 

The last decade has witnessed a tremendous progress 
in silicon-on-insulator (SOI) technology for applications 
in photonic integrated computing and data processing 
P]. In parallel, integrated photonic circuits have become 
increasingly appealing to realize key tasks in quantum 
information and communication, thanks to their natu¬ 
ral interfacing with long distance communication net¬ 
works working in telecommunication band (1.3 — 1.5 /am 
wavelengths). Clearly, the combination of these two 
paradigms will likely allow to realize complex quantum 
operations on-chip that are far beyond what may be en¬ 
visaged with table-top experiments, with significant and 
large scale impact on efficient and secure data processing 
and transmission [2]. Within this context, the genera¬ 
tion of single photons plays a central role for the devel¬ 
opment of on-chip quantum photonic technologies [3]. In 
particular, the recent advances in silicon-based quantum 
photonics Hd would strongly benefit from integrated 
single-photon sources on the same operating chip. 

Single-photon sources on-demand can be realized with 
artificial two-level emitters, such as semiconductor quan¬ 
tum dots mm, which have increasingly improved their 
radiative efficiency over the last few years jTQHT3] . How¬ 
ever, these single photon sources are typically based on 
III-V semiconductors, they work most efficiently at cryo¬ 
genic temperatures, and integration with silicon-based 
nanophotonic circuits working in the telecommunication 
band [I] remains challenging. As a possible alternative, 
integrated single-photon sources in silicon-on-insulator 
(SOI) photonic circuits have been shown [MlfTb] . based 
on enhanced four-wave mixing induced by the silicon x ^ 
susceptibility and non-deterministic heralding. Even if 
the efficiency of such integrated sources can be improved 
by spatial multiplexing m, compactness and scalability 
remain open issues. 

An alternative route to single-photon generation re¬ 


lies on the photon blockade mechanism [18], where a 
strong third-order nonlinearity in an optical resonator 
enables a shift of the resonant frequency by more than 
its linewidth when a single photon is already present. As 
a consequence, the device can absorb a photon only after 
the previous one has been emitted. However, this mech¬ 
anism however requires a stronger optical nonlinearity 
than what is achieved in state-of-the-art SOI devices m ■ 

Here, we build on the mechanism called unconven¬ 
tional photon blockade (UBP), recently advocated as a 
promising paradigm for single-photon generation [3H5]. 
The UPB mechanism relies on quantum interference, 
and is therefore highly sensitive to an optical nonlinear¬ 
ity of small magnitude. At difference with the conven¬ 
tional blockade, it has been shown that UPB can also 
occur when the nonlinear frequency shift per photon is 
much smaller than the cavity linewidth, which is usu¬ 
ally the case also in silicon photonic crystal nanocavi¬ 
ties m noii25]. So far, such theoretical mechanism was 
only shown to work under continuous-wave (cw) excita¬ 
tion, which severely limits the usefulness of the proposed 
antibunched radiation as an actual single-photon source 
HE]. In the present paper we go beyond previous works 
on UPB by demonstrating a reliable protocol that al¬ 
lows to operate any such system under pulsed excitation. 
In fact, this can be technically achieved by a combina¬ 
tion of excitation pulse tailoring and post-selective tem¬ 
poral filtering of the output stream to purify the statistics 
of the emitted radiation, similarly to what has been al¬ 
ready demonstrated for quantum dot-based single photon 
sources |T3]. The latter achievement allows to overcome 
a previously believed limitation, and provides a scheme 
to devise a true single photon source out of a generator 
of antibunched radiation. The efficiency of such a single 
photon source in realistic devices is analyzed in detail, 
which is shown to outperform the best heralded sources 
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Figure 1: Unconventional photon blockade. (a) 

Schematic representation of an asymmetrically driven pho¬ 
tonic molecule. Each cavity is characterized by a single res¬ 
onant mode in the spectral region of interest, and only the 
first cavity is driven by an external coherent field, (b) The 
corresponding ladder scheme of the lowest few energy levels, 
associated to photon occupation number states in the two 
cavities, (c) The computed second-order correlation function 
of the quantum field in the first cavity, g^ 2 \r), plotted as a 
function of time delay. The quantity was computed under the 
assumption of a cw driving field. Inset: a detail of the photon 
antibunching region close to zero delay, (d) Average photon 
occupation of the first cavity (m, full line), and corresponding 
value of g^ 2 \0) (dashed line), computed as a function of the 
cw driving field amplitude. Inset: Same quantities plotted for 
the second cavity. 


demonstrated so far in key figures of merit, especially 
operation power and spatial footprint. 

Results 

Following Refs. sn we consider UPB in a system of 
two tunnel-coupled cavities, i.e. a photonic molecule , as 
sketched in Fig. EK The quantum model of UPB has 
been thoroughly characterized [5], and it is briefly sum¬ 
marized in the Methods section. The relevant physical 
parameters are the tunnel coupling rate between the two 
cavities ( J/h ), the driving rate on the first cavity ( F/h ), 
the driving frequency (wl), the effective photon-photon 
interaction energy in each cavity (Uj, j = 1,2), which 
originates from the intrinsic material [EL an d the 
cavities loss rates Kj. Detrimental pure-dephasing pro¬ 
cesses are known to be negligible if the overall dephasing 
rate is much smaller than Uj/h\ 5]. For a photonic crystal 
cavity in silicon, this condition is largely fulfilled [261. Fi¬ 
nally, the model can be generalized to include input and 
output quantum channels [2]. We will assume ujj = cj c , 
Uj = U and Kj = k in the following, without loss of 


generality. 

A scheme of the lowest 6 levels on the basis of photon- 
number states, |niri 2 ), is given in Fig. [lj). The differ¬ 
ent excitation pathways leading from the initial ground 
state to the state 12 0) - corresponding to two-photon 
occupation of the first (driven) cavity - are highlighted. 
The UPB is essentially based on suppression of such dou¬ 
ble occupation by a careful tuning of the model parame¬ 
ters, leading to destructive quantum interference between 
the two alternative pathways. The optimal UPB condi¬ 
tions [4] are given by J op t/ha ~ [(2/3\/3)Hk/U] 1 / 2 , and 
A opt = (co c — ujl) = —k/ 2\/3, and will be assumed to 
hold in what follows. 

We consider the UPB mechanism in a SOI nanopho¬ 
tonic platform, where the cavity-field confinement in a 
diffraction-limited mode volume, V ~ (A/n) 3 , enhances 
the effective photon-photon interaction, U. A realistic or¬ 
der of magnitude estimate in a crystalline silicon photonic 
crystal nanocavity leads to U ~ lO -3 fieV (see also Sup¬ 
plementary Information) [5]. Assuming a quality factor 
Q ~ 8 x 10 5 - now routinely achieved at telecom wave¬ 
lengths (i.e., huj c ~ 0.8 eV) [6l UOl [29] - we set k ^ 1 
/ieV, and hence U/hn = 0.001. To fulfill the optimal 
UPB conditions, the remaining parameters take values 
A = —0.29 k and J = 19.6fFc, respectively. 

The steady state results under cw driving are sum¬ 
marized in Fig. [lj^-f. The time-dependent normalized 
second-order correlation function, g^ 2 \r) (see Methods 
section) is considered as the reference figure of merit for 
single-photon blockade I30H32] and plotted in Fig. [TJ^. A 
strong antibunching is present over a time-delay window 
r < 100 ps. At longer delays, strong oscillations are 
present on the timescale h/J, arising from the interfer¬ 
ential nature of the UPB mechanism [4]. The average 
photon occupations in the two cavities, rii = (aja^), and 
the corresponding zero-delay correlation, #^(0), are dis¬ 
played as a function of the driving field amplitude, F, in 
Fig. [Tp UPB occurs at low average occupation of the 
driven cavity, while the occupation of the non-detected 
cavity is much larger (see inset). This figure of merit is 
relevant to determine the maximal single-photon emis¬ 
sion rate that can be achieved in such a device under 
cw pumping, given by R em = nift/27r. As an exam¬ 
ple, for n\ ~ 0.05 (corresponding to F/Hk ~ 30 and 
g( 2 )(0) < 0.5), Rem > 10 MHz can be expected, with an 
input power as low as P[ n = hu c F/h ~ 1 nW. In fact, the 
optimal UPB relation between U and J leads to a con¬ 
dition (without numerical pre-factors, for convenience) 
J opt oc y/V/Q 2 ; this means that the required input-power 
in UPB scales down roughly as 1/Q 2 , i.e. the larger the 
cavity Q, the smaller J opt can be to have antibunching by 
keeping the average number of photons in the first cavity 
less than 0.1 (according to Fig. HP The same figure of 
merit simultaneously allows to increase the antibunching 
time window, scaling as l/J op t (see Fig- EF 

Single-photon sources on demand require the emission 
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Figure 2: Unconventional photon blockade under 

pulsed excitation, (a) Pulse sequence driving the first cav¬ 
ity, as schematically represented in (b). Pulse duration is set 
to at = 4 ns at 50 MHz repetition rate, (c) Correspond¬ 
ing cavity output, i.e. average population in the driven cav¬ 
ity as a function of time. The filtering window is schemati¬ 
cally superimposed within each pulse (yellow areas), (d) De¬ 
tail of a single pulse from the output sequence in panel (c). 
The blue curve shows the equal time second-order correlation 
function, g^ 2 \t,t) (scale on the right), (e) Un-normalized 
and un-hltered second-order correlation over the whole pulse 
sequence (see Supplementary Information), (f) Same, after 
time-filtering the pulses with AT = 75 ps. (g) Full Montecarlo 
wave function simulation of the time-filtered second-order cor¬ 
relation function under pulsed excitation, as a function of the 
width AT of the filtering window (disks), and relative error 
span (crosses). The cyan line is the two-time second-order 
correlation calculated by solving the quantum master equa¬ 
tion (see Supplementary Information). The model parameters 
assumed in these simulations are the same as in Fig. [l] 


of single-photon pulses at deterministic times. How¬ 
ever, in the UPB mechanism the emitted light is sub- 
Poissonian only within a time delay shorter than h/J , 
as shown in Fig. [ljz. For short input pulses, the outgo¬ 
ing pulses will last at least as long as the cavity lifetime. 
Therefore, the condition J Hk would apparently pre¬ 
vent the device from operating in a pulsed regime [4]. 

Here we show for the first time that UPB under pulsed 
excitation is possible by exploiting temporal filtering of 
the output signal. In Fig. [2j the results of a numerical 
experiment are reported for the UPB device considered 
in the previous section, where a train of gaussian pulses 
drives the first cavity (Fig. i'- c). A sequence of outgoing 
pulses from cavity 1 is modeled by solving the quantum 
master equation (see Supplementary Information for de¬ 
tails), and shown in Fig. Hf Focusing on a single out¬ 
going pulse, the equal-time second-order correlation is 


plotted in Fig. W (blue curve), where a well-defined time 
window clearly exists - within the pulse emitted from a 
UPB device - during which light is antibunched over a 
time delay shorter than h/J. 

As a consequence, pulsed operation can be achieved 
by gating the outgoing pulses in time, in order to retain 
only a timeframe of duration < h/J. In practice, this 
could be achieved with an integrated all-optical switch 
triggered by the input pulse, as it was already shown ex¬ 
perimentally [13]- The second-order correlation function 
under pulsed excitation (see Supplementary Information) 
is shown in Figs. !' -g. The histograms in Figs. S' -f show 
the un-normalized correlation signal, G^ 2 ^(r), integrated 
over the whole pulse sequence in (Fig. [2^) and in the pres¬ 
ence of filtering with a time window AT = 75 ps (Fig.[2]f), 
respectively. They reveal a strong reduction of the two 
photon counts within a pulse after filtering, which is a 
key result of this paper. Figure [2^ displays the depen¬ 
dence of the filtered second-order correlation versus the 
filtering time window, AT. The Montecarlo data (blue 
disks), directly obtained from the photon count statistics 
(see Supplementary Information), are reproduced by a 
master equation treatment (cyan curve), which confirms 
the reliability of this result. Photon antibunching (gray 
area) is achieved below AT = 130 ps, while the single 
photon regime - requiring the condition g( 2 \AT) < 0.5 
- is obtained for AT < 90 ps. When assuming 5 x 10 7 
pulses per second and a peak value F ~ 150/ift, after 
the temporal filtering the Montecarlo data indicate a 
single-photon yield at a rate of about 0.45 MHz. Under 
these conditions, the driven cavity reaches a peak value 
of the average photon occupation n\ ~ 0.075, close to the 
largest occupancy for which UPB is expected according 
to Fig. ^ [i.e., g( 2 )(0) < 0.5]. Remarkably, this peak 
value of F(t) implies an intracavity energy of less than 
10 -2 fJ per pulse. This corresponds to an input energy 
that can be quantified as ~ 0.5 fJ per pulse, according to 
typical excitation schemes of photonic crystal integrated 
circuits [25] . We notice that this is extremely low when 
compared to the state-of-art parametric sources demon¬ 
strated so far in integrated silicon-based platforms and 
based on four-wave mixing and heralding (typical input 
powers in the 100 mW range), for a comparable output 
rate in the few 100 kHz range mm- 

A feasible realization is hereby proposed in an in¬ 
tegrated SOI photonic crystal platform. These struc¬ 
tures benefit from a remarkably advanced and well es¬ 
tablished fabrication technology, with nanocavities hav¬ 
ing recently achieved Q/U-values well above the UPB 
requirements 0 EEl E3 ESt As a schematic example, 
we show in Fig. [3^i a double-cavity device in a photonic 
crystal circuit. This configuration allows to selectively 
drive one cavity from the input waveguide channel, and 
simultaneously collect part of the light emitted from the 
same cavity into the output waveguide (the remaining 
part being emitted through out-of-plane losses). As an 
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Figure 3: Realization of a SOI integrated single¬ 

photon source, (a) Artistic view of an integrated SOI pho¬ 
tonic crystal chip realizing input/output channels and UPB 
through a photonic crystal molecule, (b) Fine-tuning of the 
normal mode splitting of a L3 photonic crystal molecule (see 
text) through variation of the radius (7*2) in the middle row 
of holes (red-highlighted in the inset). The holes highlighted 
in green are shifted off from the cavities center to optimize 
the Q-factor ehid]. The hole radius of the surrounding pho¬ 
tonic crystal lattice is r — 112 nm. All the design parameters 
are given in detail in the Supplementary Information associ¬ 
ated to this manuscript, (c) Spectrum of the photonic crystal 
molecule designed to have the parameters corresponding to 
the results shown in Figs. [Tp] The two normal mode reso¬ 
nances are indicated as bonding (B) and antibonding (AB), 
respectively. The exciting laser frequency at the optimal UPB 
condition is schematically indicated (vertical dashed line). In¬ 
set: E y component plotted for the two normal modes, super¬ 
imposed to the photonic crystal design showing the footprint 
area of the coupled cavities device. The reference directions 
x and y are explicitly indicated. 


elementary building block, we consider a L3 photonic 
crystal cavity in a thin silicon membrane, designed for 
operation at the preferred telecom wavelength, A = 1.5 
/im 0.825 eV). The cavity consists of three missing air 
holes in a triangular air-hole lattice. This cavity was re¬ 
cently optimized to show a measured quality factor reg¬ 
ularly exceeding one million UM- We use here a L3 
cavity design with theoretical unloaded (i.e., valid for 
the isolated cavity) Q ^ 1.25 x 10 6 (see Supplementary 
Information for details on the structure parameters, such 
as hole radius and lattice constant). When coupling to 
the access waveguides, the loaded Q-factor can be engi¬ 
neered in the range Q ~ 8 x 10 5 , as we have verified by 


3D finite-difference time-domain simulations (3D-FDTD, 
not shown). From the calculated mode profile, the effec¬ 
tive nonlinearity for this device is estimated (see Sup¬ 
plementary Information) in the range U ~ 0.8 x 10 -3 
//eV, close to what was assumed in the model calcula¬ 
tions above. 

The photonic crystal molecule can be obtained by ver¬ 
tically aligning two L3 cavities, separated by 5 rows of 
holes (i.e., 3^/3a center-to-center) [33]. The hole radius 
in the central row, 7 * 2 , can be varied to fine tune the nor¬ 
mal modes splitting at the desired value [53] * i.e. A = 2 J 
according to Eq. [I] in Methods. In Fig. [ 3 J 3 we show such 
a simulated fine tuning by 3D-FDTD calculations. For 
7*2 ^ 95 nm, the normal mode splitting between the two 
resonances, identified as bonding (B) and antibonding 
(AB) according the the spatial profile of the E y compo¬ 
nent, is calculated as A ~ 44 /ieV, i.e. remarkably close 
to the condition J/hn = 19.6 assumed in the previous 
calculations, when we consider the loaded value Kk ~ 1 
/ieV. We notice that similar values and dynamic control 
of the normal mode splitting have been already shown ex¬ 
perimentally in SOI photonic crystal platforms operating 
in a very similar wavelength range [35]. The spectrum 
for such an optimal structure is shown in Fig. [ 3 J 3 . The 
two resonances have unbalanced Q-factors of ~ 1.1 x 10 6 
(AB) and ~ 1.3 x 10 6 (B), respectively, which can also be 
exploited to enhance the degree of antibunching in UPB 
0 . 

Before concluding, we discuss how to circumvent the 
most relevant and potentially detrimental effects for the 
realization of UPB in a SOI platform. First of all, two- 
photon absorption (TPA), related to the imaginary part 
of the silicon is also enhanced by confinement in 

the L3 cavities. However, a quantitative estimate of 
this contribution has been given in Ref. [13] by which 
we can infer a TPA loss rate that is on the order of 
ktpa/k < 10 -2 for the present case, also considering the 
low input powers necessary to achieve UPB. Second, ther¬ 
mal effects can give rise to pure dephasing of the cavity 
resonances, which depends on optomechanical coupling 
with the background phonons. For the L3-type silicon 
photonic crystal cavities considered here, this contribu¬ 
tion has been quantitatively estimated and shown to be 
negligible even at room temperature (i.e., a pure dephas¬ 
ing rate 7*/ft ^ 10 -7 ) [26] . Finally, unavoidable fabri¬ 
cation imperfections should be corrected by device post¬ 
fabrication processing. In particular, fine and selective 
cavity tuning has been already shown for photonic crys¬ 
tal cavities with different techniques, even in the presence 
of very large Q-factors [35H37] . 

Discussion 

We have theoretically shown that an integrated 
nanophotonic platform based on CMOS-compatible SOI 
technology can be engineered to achieve single photon 
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emission by an unconventional photon blockade mecha¬ 
nism. Besides opening the way to the first experimental 
demonstration of UPB, our results show that such uncon¬ 
ventional mechanism allows for pulsed excitation, which 
represents a key ingredient for a useful source where each 
pulse potentially triggers emission of a single photon. 

Such an alternative single-photon source could be char¬ 
acterized by a very low input power operation, i.e. com¬ 
parable to standard single-photon devices based on cav¬ 
ity QED but much lower than typical integrated single¬ 
photon sources based on enhanced four-wave mixing and 
heralding. Moreover, this is achieved by an unprece¬ 
dented small footprint area, significantly smaller than re¬ 
cently realized heralded sources in integrated SOI chips. 
In fact, notice that the footprint of this prospected device 
is essentially given by the spatial extension of the pho¬ 
tonic crystal molecule and the necessary lattice around it. 
For the structure simulated in Fig. |3j we estimate a mini¬ 
mal footprint area on the order of a few fim 2 (see, e.g., the 
inset in Fig. [3^). In practice, this is significantly smaller 
than current heralded sources fabricated with the same 
SOI technology and based on coupled resonator optical 
waveguides M or spatially multiplexed photonic crystal 
waveguides D2I- We also stress the generality of the pro¬ 
posed scheme, which could be extended to other types of 
nonlinearities [38], and could eventually lead to the real¬ 
ization of novel quantum devices (13 HO] for applications 
in integrated quantum metrology and logic. 

In summary, by combining an extremely low input 
power, a small footprint area, and no quantum emitter re¬ 
quired for single-photon generation, such a device might 
have significant impact on the development of integrated 
silicon quantum photonics, by introducing a new concept 
in the generation of pure quantum states of light at ar¬ 
bitrary wavelengths (e.g., in the telecom band), that is 
fully compatible with current semiconductor technology, 
working at room temperature, and a viable alternative 
to single-photon nonlinear devices based on cavity-QED 
with artificial atoms or single atomic-like emitters that 
are presently lacking in SOI integrated platforms. 

Methods 

Model Hamiltonian. The second quantized hamiltonian of the 
driven nonlinear photonic molecule is expressed (to leading 
linear and nonlinear orders) as [3] 5. 

7 is = ^ [hjjjd\di + Uid\d\djdi ] + J(d\d 2 + a\di) 

£= 1,2 

+ F(t)e~ iULt a{ +F*(t)e i “ Lt a 1 . (1) 

The first terms in Eq. ^ describe two harmonic oscillators, 
J/h is the tunnel coupling rate between the two resonators, 
F(t)/h is the coherent pumping rate on the first cavity at the 
laser frequency c ol, and the photon-photon interaction en¬ 
ergy in each resonator is related to the material x^ min]- 
A description of this quantity and an estimation for the pho¬ 
tonic crystal cavities considered here are given in the Sup¬ 
plementary Information. Starting from this Hamiltonian, the 


various time-dependent photon correlation functions for light 
collected after cavity 1, generally defined as 

(2) /, ./>. (a!(t)al(Qai(t')ai(t)) 

’ (oI(t)o 1 (*))(a t 1 (*')o 1 (t')) ’ 

were numerically simulated by using both the Montecarlo 
wave function method and by directly solving the master 
equation for the density matrix (see details in Supplemen¬ 
tary Information). In both cases, the numerical solutions 
were computed on a truncated Hilbert space of dimensions 
(N -_l x N 2 ) 2 , where N± = 4 and N 2 = 18 are the maximum 
photon occupations allowed in cavities 1 and 2, respectively. 
While the master equation results were obtained from a mod¬ 
ern workstation embedding 16 Gb of RAM, the Montecarlo 
data were produced by 10 nodes of 16 cores and 32 Gb RAM 
memory each, run on a high-end cluster for a few weeks of 
continuous computational time. 
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We detail here the Montecarlo and master equation treatments used to produce the results of figures 1 and 2 in 
the main text. We also give details on the photonic crystal cavities design, and the corresponding estimation of the 
single-photon nonlinearity reported in the manuscript. 


MONTECARLO WAVE FUNCTION METHOD 

The statistics of the photons emitted by the system 
under pulsed excitation were first addressed by perform¬ 
ing quantum Montecarlo simulations [l|. This method 
not only allows to work with larger truncated Hilbert 
spaces but also provides direct access to individual pho¬ 
ton counts, thus embodying the closest theoretical sim¬ 
ulation of an actual Hanbury Brown-Twiss experiment. 
In brief, the algorithm is based on the stochastic evolu¬ 
tion of the system wave function through the Schrodinger 
equation 

^ i^> * ih §t i^> ’ w 

written for the non-Hermitian effective Hamiltonian 

l~i = 1~L S — — ^ tvjOjdj , (2) 

3 

The non Hermitian part of [2] results in a decay of 
the norm | i/). During the evolution of Eq. [I] ran¬ 
dom numbers 0 < r < 1 are drawn and the condition 
(-0 | 'ip) < r decides for the action of a jump operator, 
dj |^), corresponding to the measurement of a photon. 
The proper quantum jump operator is chosen such that 
j is the smallest integer satisfying JT Pj > r, where Pj 
are the probabilities for the mode j to emit a photon at a 
given time. Each evolution of Eq. [l] produces a stochas¬ 
tic quantum trajectory associated with the state \ip(t))j, 
and the procedure can be repeated N times to form an 
ensemble average of realizations in view of approximat¬ 
ing the system density matrix as p(t) = 14/ (£))(4/ (£)I, 

N-t oo 

where 

N 

\*(t)) = '£\^(t)) j /N, (3) 

3 =1 

and its potential mixed nature. Any observable or cor¬ 
relation are obtained from ( 0(t)) = Tr [p(t)0]. The full 
access to photon counts and emission times history al¬ 
lows to mimic the experimental detection scheme. The 
two-times second-order correlation function, #( 2 )(£i,£ 2 ), 
can be reconstructed from the statistics of photons de¬ 
lays analogously to a Hanbury Brown and Twiss (HBT) 
setup. Further details on the numerical procedure em¬ 
ployed to obtain the results of Fig. 2 in the main text are 
given in the following. 


TWO-TIME CORRELATIONS UNDER PULSED 
EXCITATION 

In our Montecarlo simulations we have worked on the 
basis of trajectories containing single pulses, which makes 
the data analysis more flexible. We have tracked the 
quantum jumps performed by the driven cavity from 
which the photon antibunching is expected. We point out 
that in reality one should expect a weak mixing between 
both cavity fields to occur in the guiding channels. It can 
be accounted for within an input-output treatment. In 
such a case the laser detuning should be properly adapted 
according to the prescriptions of Ref. SJ To analyze de¬ 
lays within a given pulse we need to track the trajectories 
where at least two quantum jumps occurred within AT 
and these events are obviously rare given the relatively 
small occupation of the cavity 1 as one can see from Fig. 2 
of the main text. We have therefore performed a large 
campaign of massively parallelized simulations on an high 
end cluster based on N = 1.5 x 10 8 pulses from which 
we have recorded the whole quantum jump history. We 
considered pulses of duration 4 ns separated by 20 ns 
to avoid any overlap bringing some unwanted pulse to 
pulse correlations. Our simulation therefore covers not 
less than 1 seconds of recorded data. 

To build the Monte-Carlo curve of Fig. 2g (blue disks), 
we have worked on quantum jumps that occurred in a 
time window of width AT = T 2 — T\ = 1.57 ns centered 
on the g( 2 \t,t) minimum (see yellow surface and blue 
curve of Fig. 2f) mimicking the temporal filtering. In¬ 
side this global window we have considered sub-windows 
of variable duration At = £2 — £1 ranging from 6.5 ps to 
AT. Each of these sub-windows was gradually displaced 
by At within AT, starting from the condition t\ = T\ 
(£2 = £1 + At) and until £2 = T 2 (£1 = £2 — At) is fulfilled. 
For a given value of A£, the un-normalized second order 
correlation G^ (At) is obtained from the sum of photon 
pair counts recorded by slicing the time window, which 
increases the statistics by N w = AT/At. Therefore it 
allows to work with a number of counts that would cor¬ 
respond to Ax N w trajectories (pulses) reducing by N w 
the required computational time. Finally, the g( 2 \At) is 
obtained by normalizing G^ (At) to the number double 
counts expected from a Poissonian statistics. The errors 
(magenta curve) are computed from the square root of 
the number of counts averaged over the sliding windows. 
Obviously the error is inhomogeneous versus A£, given 
that N w is variable, and it is small both in the regions of 



narrow and wide At where respectively N w is large and 
number of counts is important. The previously described 
procedure is summarized in Fig. SI (see captions). 

To build the histograms of Fig. 2d-e displaying the 
pulse-to-pulse statistics, we have performed a Monte- 
carlo rearrangement of our single pulse trajectories to 
randomize their time ordering, as it would be obtained 
from many pulse trajectories or in an actual experimental 
situation. 



Figure 1: (a) Average occupation of cavity 1 (red) and the 
equal time second order correlation function (blue line). The 
yellow region highlights the global time window AT from 
which the photon counts are extracted, (b) Zoom in between 
the vertical dashed lines of panel (a) showing an illustrative 
set of sub-windows for a given value of At. The wavy lines 
illustrate and example of quantum jump series for a given 
trajectory. Only the 2 jumps event are kept for the G^ 2 \At) 
statistics. 3 and more jumps events are totally absent in the 
conditions we consider. 


QUANTUM MASTER EQUATION 

The master equation for the density matrix reads 

p = + £ {1,2) , (4) 

where losses are accounted for through Liouvillian oper¬ 
ators in the usual Lindblad form for the two resonators 
modes, £ (1 ’ 2) = Ki[aipa\ - 0.5(aja*p - pat a*)]. 

Further sources of loss, such as nonlinear absorption (e.g. 
related to the imaginary part of x^) or pure dephasing, 
could also be added to Eq. [4] (see, e.g. [5j [13]), but we 
neglect them here for simplicity. Moreover, the model 
can be generalized to include input and output quantum 
channels [2|. 

Single-time evolution and steady state numerical re¬ 
sults of Eq. [4] can be straightforwardly performed, as 
in Refs. IMl Here, we were additionally able to con¬ 
firm the Montecarlo results (cyan curve in Fig. 2g of 
the main text). The un-normalized and normalized two- 
times second-order correlation functions of cavity 1 were 


computed as 


g (2) (m') 
g {2) (t,0 


Tr 


a\a±U t ^t' (al«i/3(0) 


2) (M') 


Tr 


a\aip(t) Tr a\a\p{t') 


( 5 ) 

(6) 


where is the propagator of the operator O from 

t to t' associated with Eq. [4] The photon statistics pro¬ 
duced within a time window At = t 2 — t\ is obtained 
from 


(2) (At) = //aT (2) (t,f)dtdf 
I /a* Ul (t')dtdt' 


( 7 ) 


where n\(t) = Tr[a{dip (£)]. This exact calculation per¬ 
fectly reproduces the Montecarlo wave function results 
within the error envelope, as it is reported in Fig. 2g 
(cyan curve). 


PHOTONIC CRYSTAL CAVITIES DESIGN 


Photonic crystal cavities allow to achieve record fig¬ 
ures of merit today, such as ultra-small mode volumes 
and ultra-high quality factors [6]. One of the most used 
photonic crystal cavity designs is realized by removing 
three air holes in a triangular lattice [T, which is usu¬ 
ally defined a L3 point defect. Recently, a combination 
of fast simulation tools [5j and genetic optimization |9j 
have allowed to show that Q factors largely exceeding 
10 6 can be designed for such cavities, which was shown 
experimentally m- 

For the photonic crystal cavities design used in this 
work, we started from a standard SOI photonic crystal 
membrane, with the silicon layer thickness of 220 nm. 
We set the lattice constant to a = 400 nm and the holes 
radius to r = 112 nm ( r/a = 0.28) to tune the cavity 
mode resonant wavelength in the relevant telecom band, 
i.e. A = 1.5 fi m (^ 0.825 eV). The three holes along the 
cavity axis have been shifted by si = 120 nm ( 51 /a = 0.3) 
5 2 = 100 nm (s 2 /a = 0.25), and 53 = 40 nm ( 53 /a = 0.1), 
to reach a theoretical (unloaded) Q 1.25 x 10 6 . Since 
we aim at coupling these cavities with access and output 
waveguides, we thus allow the loaded Q-factor to be in 
the 0.8 x 10 6 range. 

A plot of the normalized electric field intensity, i.e. 
.the function |^(r)| 2 with J|<a(r)| 2 dr = 1, is shown in 
Fig. S2 for our optimized L3 cavity design, which was 
the building block for the the photonic crystal molecules 
in Fig. 3 of the main text. 
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Figure 2: Electric field intensity profile at the center of the 
silicon photonic crystal membrane. 


ESTIMATING THE EFFECTIVE 
PHOTON-PHOTON INTERACTION 


In this work, we focus on photonic nanostructures in 
silicon, which is a strongly nonlinear material already at 
the level of classical electromagnetic response. In par¬ 
ticular, bulk silicon is characterized by a relatively large 
X< 3) susceptibility, while nominally x^ = 0 (neglect¬ 
ing surface contributions) owing to the centrosymmetric 
nature of the elementary crystalline cell m- Strongly 
enhanced nonlinear effects have been already reported in 
L3 photonic crystal cavities [12 . 

The photon-photon interaction energy in each res¬ 
onator is given in terms of the material x^ by the sim¬ 
plified expression m 


u = 


D(hu)i) 2 f x (3) ( r ) 
8e 0 J r e 2 (r) 


Kr)| 4 , 


(8) 


where a(r) is the three-dimensional cavity field profile, 
normalized as f |cf(r)| 2 dr = 1, and D represents the mul¬ 
tiple contributions of the same order of magnitude given 
by the different elements of the x® tensor m- 

From the calculated mode profile shown in Fig. S2, the 
effective nonlinearity of such a silicon nanocavity can be 
estimated through Eq. [5J by using x® ^ 0-9 x 10“ 18 
m 2 /V 2 , which is an appropriate order of magnitude for 
the elements of the bulk silicon third-order susceptibility 
tensor [Tl], and D=24 [14]. For the cavity mode profile 
of our optimized photonics crystal structure, see Fig. S2, 
a quantitative estimate of this integral results in U ~ 
0.8 x 10 -3 /ieV, close to what was assumed in the model 


calculations of the main text and confirming the order of 
magnitude estimates already given in the literature [5]. 
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